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ABSTRACT 

This  paper  describes  and  analyzes  a  Gaussian  statistical  satellite  detection  algorithm.  This  detection  algorithm 
exploits  the  correlation  of  pixels  that  arises  from  the  long-term  exposure  point  spread  function  (PSF)  due  to 
diffraction  and  atmospheric  effects.  This  correlating  detection  algorithm  is  compared  to  a  simple,  linear  threshold 
detection  algorithm,  which  treats  the  probability  density  function  (PDF)  of  the  photons  as  Gaussian  and  the  PSF  as 
Dirac.  The  effects  of  changing  the  window  size  and  grid  size  is  explored  and  analyzed.  Using  a  dataset  collected  by 
the  Space  Surveillance  Telescope  (SST)  that  shows  many  samples  of  a  geostationary  satellite  gradually  going  into 
eclipse  during  the  vernal  equinox,  the  probabilities  of  detection  for  the  algorithms  are  compared  as  the  satellite 
becomes  very  dim  and  nearly  disappears  to  determine  which  detection  algorithm  performs  better.  A  better 
performing  detection  algorithm  will  allow  detection  of  satellites,  space  debris,  and  other  dim  objects  that  have  not 
been  visible  before. 


1.  INTRODUCTION 

This  paper  investigates  what  gains  can  be  achieved  by  using  a  non-linear  log-likelihood  detector  that  exploits  the 
correlation  that  results  from  diffraction  over  current  techniques  that  only  consider  the  light  from  a  single  pixel.  It 
also  analyzes  the  effects  of  changing  the  size  of  both  the  window  and  the  search  grid  on  the  performance  of  both 
algorithms.  First  the  experiment  that  provided  data  for  this  analysis  is  introduced.  Then  a  mathematical  Gaussian 
detection  algorithm  is  provided,  followed  by  a  discussion  of  the  methodology  and  limitations  of  the  comparison. 
Finally  the  results  are  discussed  and  analyzed.  That  correlation  provides  a  better  probability  of  detection  is 
expected.  Without  an  understanding  of  the  conditions  under  which  both  the  single  point  and  the  correlated  detectors 
work  best,  though,  it  is  difficult  to  determine  whether  either  is  being  implemented  wisely  in  a  given  situation. 

2.  EXPERIMENT 

The  goal  of  this  research  project  is  to  compare  the  performance  of  detection  algorithms  under  very  low  light 
conditions.  A  glint  experiment  done  by  the  USNO  (United  States  Naval  Observatory)  this  past  fall  served  as  the 
blueprint  for  the  experiment  that  we  conducted  [1].  In  their  experiment,  a  geostationary  satellite  was  imaged  going 
into  eclipse  during  the  autumnal  equinox.  Data  from  the  USNO  in  table  1  shows  how  a  rapid  series  of  images 
revealed  the  satellite  growing  dimmer  gradually  until  it  disappears  entirely.  This  situation  presents  the  observer  with 
an  object  whose  intensity  is  varying  in  time,  which  allows  detection  algorithms  to  be  tested  over  a  variety  of  light 
levels.  Also  because  the  object  can  be  clearly  observed  heading  into  the  eclipse,  its  position  is  relatively  well 
known  throughout  the  experiment.  These  attributes  make  this  type  of  test  ideal  for  comparing  the  performance  of 
different  detection  algorithms. 
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Table  1.  Digital  counts  for  Anik  FI  entering  the  terminator  during  autumnal  equinox.  Provided  by  Dave  Monet, 
USNO 


From  USNO  Telescope 
w/Orthogonal  Transfer  CCD* 

Frame 

UT 

Counts  above  sky 

2500 

25:31.1 

saturation 

2550 

40:41.1 

7100 

2560 

41:31.0 

2000 

2570 

42:21.2 

1100 

2580 

43:11.2 

700 

2590 

44:01.3 

400 

2600 

44:51.2 

200 

2610 

45:41.2 

200 

2620 

46:31.4 

200 

2630 

47:21.4 

150 

2640 

48:11.2 

50 

2650 

49:01.3 

0 

In  our  experiment  the  recently-constructed  SST  telescope  imaged  a  geostationary  communications  satellite  as  it 
went  in  and  out  of  eclipse  for  eleven  nights  between  28  February  and  23  March  2012  during  the  vernal  equinox.  We 
also  imaged  Anik  FI,  which  flies  at  107.3°  W  longitude  and  is  40.4  x  9  meters  at  its  largest  dimensions.  The 
satellite  was  actively  tracked  from  its  orbital  elements  through  the  entire  eclipse.  The  result  is  a  dataset  in  which  a 
very  dim  object  is  in  a  known  location,  providing  the  opportunity  to  test  detection  algorithms  under  a  wide  variety  of 
light  conditions.  Since  the  satellite  was  tracked,  other  objects  in  the  sky  pass  through  at  the  sidereal  rate.  Pictures 
were  taken  every  five  seconds,  resulting  in  relatively  gradual  changes  between  frames.  Both  binned  (2x2)  and  un¬ 
binned  pixel  data  was  collected,  but  only  binned  data  is  included  in  this  analysis. 

The  dataset  does  not  provide  a  basis  for  forming  receiver  operating  characteristic  (ROC)  curves.  The  design  of  SST 
is  unprecedented,  and  it  is  certainly  capable  of  seeing  things  that  have  not  been  seen  before.  Thus,  it  is  extremely 
difficult  to  assert  that  it  is  imaging  a  blank  sky,  which  is  the  necessary  condition  for  computing  the  probability  of 
false  alarm. 


3.  GAUSSIAN  LOG-LIKELIHOOD  DETECTORS 

3.1  Mathematical  Basis 

The  satellite  is  considered  to  be  present  in  a  pixel  d  (x,  y )  in  the  image  when  the  likelihood  ratio,  A,  exceeds  a 
threshold,  77  ,  as  in: 


A(x,  y) 
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P[d{w,z)\H0]< 


The  model  that  assumes  that  the  noise  in  the  image  is  Gaussian  is  considered  first.  The  probabilities  that  the  object 
is  present  (  //| )  or  not  (  H  )  are: 
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where  d(w,z )  are  pixels  in  a  square  NxN  window  of  the  image  centered  on  d(x,  y) ,  B  is  the  expected  value  of 
the  background  of  the  window,  <7  is  the  standard  deviation  of  the  background  of  the  window,  6  is  a  scale  factor  for 
the  expected  amount  of  light,  and  li  is  the  point  spread  function  of  the  system.  Whereas  analysis  of  data  from  other 
experiments  like  the  Panoramic  Survey  Telescope  &  Rapid  Response  System  (Pan-Starrs)  is  able  to  make  use  of 
foreknowledge  of  the  background,  the  background  is  approximated  in  this  comparison  by  the  mean  of  the  image  in 
the  window: 
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This  is  a  good  approximation  when  the  window  size  is  chosen  to  include  the  tails  of  the  PSF  (point  spread  function) 
so  that  the  intensity  of  most  of  the  pixels  in  the  window  is  dominated  by  the  background. 

Likewise,  er  ,  the  standard  deviation  of  the  background,  is  given  by: 
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Taking  the  natural  log  changes  the  products  to  summations,  allows  the  removal  of  constants  and  a  lower  dynamic 
range  of  equation  [2],  Simplifying,  the  result  is: 
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Since  d(w,  z)  is  assumed  to  be  Gaussian,  then  when  the  background  is  subtracted  and  that  sum  is  divided  by  the 
standard  deviation,  the  result  is  a  normal  Gaussian  variable,  with  zero  mean  and  unit  variance.  [1] 


3.2  Transfer  Functions 

Examining  only  the  center  pixel  with  coordinates  (cx ,  c  )  and  ignoring  correlation  in  the  rest  of  the  image  is 

equivalent  to  treating  the  transfer  function  of  the  system  as  a  Dirac  delta  function.  Accordingly,  we  can  use  the 
sifting  property  of  the  Dirac  to  make  the  following  simplification  for  the  Gaussian-Dirac  detector: 


1nA(cv,C,)  =  XZ 
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In  contrast,  correlation  in  the  image  can  be  included  by  using  the  PSF  as  the  transfer  function  of  the  system.  Thus, 
the  Gaussian-PSF  detector  considers  the  correlation  of  all  the  pixels  in  the  window: 
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lnA(cr,cy)  =  2,Z. - 1 - - 

»•  * 

4.  METHODOLOGY 


4.1  Determination  of  the  PSF 

The  normalized  PSF  of  the  satellite  for  a  given  eclipse  event  is  derived  from  dpsF  (x,  y) ,  an  average  of  the  first  ten 
frames  after  the  satellite  stops  glinting 
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where  d(x,  y)  are  images  of  the  satellite  and  (x,  y)  are  elements  of  an  /V  x  TV  window  centered  on  the  satellite.  The 
normalized  PSF  is  computed  by 


,  dPSF(x,y)-mm(dPSF) 
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Subtracting  min(<f  ) ,  which  is  the  lowest  value  found  in  dpsF  (x,  y) ,  removes  as  much  of  the  background  as 

possible  while  preserving  the  tails  of  the  PSF.  Foreknowledge  of  the  PSF  of  the  satellite  certainly  helps  in  detection. 
Essentially,  it  helps  to  know  what  it  is  that  one  is  looking  for,  and  there  are  many  cases  in  which  one  does  have  an 
opportunity  to  image  an  object  before  subsequently  searching  for  it.  A  generalized  PSF  could  also  be  used  that  was 
based  on:  the  parameters  of  the  telescope  used;  the  long  exposure  OTF  (optical  transfer  function)  of  the  atmosphere; 
and  the  general  knowledge  of  the  size  and  shape  of  the  object  being  looked  for. 

4.2  Search  Pattern  and  Logic 

The  image  of  the  satellite  moves  perceptibly  from  frame  to  frame.  Reasons  for  this  include  atmospheric  diffraction 
(mainly  tilt  aberrations),  errors  in  tracking,  orbital  path,  and  changes  in  how  the  light  is  reflected.  Therefore,  the 
algorithms  are  tried  on  all  the  pixels  in  a  square  search  grid  that  is  roughly  centered  on  the  satellite.  The  satellite  is 
considered  to  be  detected  when  it  is  found  in  any  pixel  in  the  search  grid,  the  equivalent  of  applying  a  logical  “or.” 

For  the  purposes  of  this  analysis,  the  eclipse  is  considered  to  last  from  the  first  non-saturated  frame  to  the  last  non- 
saturated  frame.  The  algorithm  was  only  used  to  detect  objects  during  this  defined  eclipse  period  for  a  couple  of 
reasons.  First,  the  objects  are  easily  detectable  in  saturated  frames  by  simply  searching  for  pixels  that  have  an 
intensity  that  equals  the  maximum  value  of  the  ADC,  so  it  is  a  trivial  case.  Second,  at  saturation  charge  spillover 
into  nearby  pixels  distorts  the  data  and  saturation  is  in  itself  a  significant  non-linearity,  both  of  which  challenge  the 
assumption  that  the  system  as  a  whole  is  essentially  linear. 

The  threshold  rj  is  set  to  six  for  both  algorithms,  so  that  they  are  competing  on  an  even  basis.  This  was  chosen 
because  LINEAR  (Licoln  Near-Earch  Asteroid  Research)  has  been  employing  a  threshold  of  77  =  6 in  their  work  [3]. 
For  the  purposes  of  comparison,  the  relative  intensity  b(x,  y)  of  a  given  point  approximates  the  intensity  of  the 
noiseless  signal  and  is  defined  as 


b(x,  y)  =  d  (x,  y)  —  B{x,  y) 


where  d(x,  y)  is  the  intensity  of  pixel  and  B(x,  y)  is  the  background  at  that  same  point,  as  defined  previously. 


5.  LIMITATIONS 


To  simplify  the  calculations,  the  windows  used  to  compute  statistics  and  to  contain  the  PSF  are  square,  even  though 
the  diffracted  image  is  better  approximated  as  a  circle.  If  the  window  size  is  chosen  well,  then  this  effect  should  be 
minimized  since  the  tails  of  the  sine  function,  which  well  describes  diffraction  [4],  fall  off  quickly.  However,  any 
square  window  will  have  too  much  noise  in  the  corners  and  not  enough  of  the  tail  along  the  side. 

If  the  window  size  is  too  large  or  small  relative  to  the  size  of  the  object  present,  then  the  statistical  detector  cannot 
distinguish  the  object  from  the  background.  This  problem  is  especially  notable  with  the  Gaussian-Dirac  detector. 
For  example,  when  using  a  33-pixel  wide  window  to  process  a  5  by  5  grid  centered  on  the  satellite  right  before  the 
beginning  of  the  eclipse  on  the  night  of  2  March,  2012  the  Gaussian-Dirac  detector  does  not  detect  the  satellite  at  all 
until  the  264th  frame,  whereas  the  Gaussian-PSF  detector  finds  the  satellite  from  the  very  first  frame;  the  first  non- 
saturated  image  is  frame  281. 


6.  RESULTS 

The  total  number  of  detections  for  the  two  competing  algorithms  on  five  nights  is  shown  in  fig.  1 .  Odd  window 
sizes  between  1  and  49  were  tried.  The  results  of  different  search  grid  sizes  are  shown  for  2  Mar  2012.  In  terms  of 
total  detections,  the  Gaussian-Dirac  never  performed  as  well  as  the  Gaussian-PSF.  Frame  by  frame  results  are 
shown  in  fig.  2  for  each  of  the  tests.  The  character  of  2  Mar  2012  is  clearly  different  from  the  other  nights,  as  the 
satellite  seems  especially  apparent  to  the  Gaussian-PSF  detector  long  after  it  disappears  from  the  Gaussian-Dirac 
detector.  Some  of  the  detections  that  are  well  after  the  eclipse  has  begun  seem  to  be  due  to  star  crossings,  but  others 
seem  legitimate;  for  most  of  those  late  detections  it  is  difficult  to  say  whether  a  dim  star  is  present  in  the  background 
or  whether  the  satellite  is  truly  detected. 

The  rapidity  with  which  the  satellite  enters  the  eclipse  changes  from  night  to  night.  In  an  attempt  to  normalize  for 
this  difference,  it  is  useful  to  examine  the  relative  intensity  b(x,  y )  of  the  data  for  detections  near  the  transition  into 
total  eclipse,  after  which  the  satellite  seems  to  be  essentially  undetectable.  The  selection  of  these  points  was  done  by 
hand  by  choosing  the  successful  detection  near  the  end  that  seemed  to  have  the  lowest  light  level.  For  the  night  of  2 
Mar  2012,  frame  500  was  the  last  frame  with  a  detection  by  the  Gaussian-Dirac  algorithm,  and  532  was  the  last 
frame  with  a  detection  by  the  Gaussian-PSF  algorithm.  These  frames  were  examined  to  determine  which  pixels 
registered  the  detection,  and  b(x,  y)  was  computed  for  those  pixels.  The  minimum  of  these  pixels  for  a  given  frame 
is  represented  by  b-psf  for  the  Gaussian-PSF  and  b-Dirac  for  the  Gaussian-Dirac.  In  computing  the  background 
B  for  use  in  finding  b(x,  y )  for  each  frame  and  algorithm,  the  window  size  that  detected  the  most  frames  in  that 
dataset  was  used  (the  maximum  value  for  each  algorithm  in  fig.  1). 

The  ratio  of  their  intensities  was  found  to  be 


bpSF  (*>  y) 


bDiraM  y ) 


0.374 


which  can  be  interpreted  to  mean  that  the  Gaussian-PSF  detector  found  the  satellite  when  it  was  only  about  a  third  as 
bright. 

For  the  night  of  14  Mar  2012,  frames  834  and  854  were  selected  for  the  Gaussian-Dirac  and  Gaussian-PSF 
detectors,  respectively;  although  834  was  not  the  last  frame,  it  was  certainly  the  frame  with  the  lowest  light  in  the 
adjacent  area.  The  ratio  of  their  intensities  was  found  to  be: 
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Fig.  1.  Comparison  of  the  number  of  detections  for  both  algorithms.  The  Gaussian-Dirac  is  shown  in  blue  in  front. 
The  Gaussian-PSF  is  shown  in  green  and  is  behind  the  Gaussian-Dirac.  The  algorithms  were  computed  using  odd 
window  sizes  between  1  and  49  on  each  of  the  different  nights.  The  search  grid  is  9x9,  unless  otherwise  noted. 
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Fig.  2.  Detections  of  the  satellite,  frame  by  frame,  for  the  window  size  that  yielded  the  most  detections  according  to  the 
total  number,  bounded  by  49x49.  If  multiple  window  sizes  were  equally  effective,  the  smallest  was  used.  The  upper  third 
shows  frames  that  only  the  Gaussian-Dirac  algorithm  detected  the  satellite.  The  middle  third  shows  frames  that  both  the 
Gaussian-Dirac  and  the  Gaussian-PSF  algorithms  detected  the  satellite.  The  lower  third  shows  frames  that  only  the 
Gaussian-PSF  algorithm  detected  the  satellite.  The  search  grid  was  9x9,  unless  otherwise  noted.  The  initial  frame  is  the 
first  non-saturated  frame. 


7.  ANALYSIS 


Past  a  certain  minimum,  the  size  of  the  window  has  a  relatively  small  impact  on  the  results,  especially  in  proportion 
to  the  additional  computational  requirements.  This  shelving  is  reached  more  quickly  and  is  more  pronounced  for  the 
Gaussian-PSF  detector  than  it  is  for  the  Gaussian-Dirac  detector.  The  best  explanation  for  the  difference  in  the 
results  for  2  Mar  2012  from  the  other  nights  seems  to  be  that  for  whatever  reason  the  satellite  seemed  to  pass  more 
slowly  into  the  terminator  on  2  Mar  2012.  Though  that  night  does  seem  to  be  a  bit  of  an  outlier,  the  smoothness  of 
its  curves  in  fig.  1  lend  credence  to  the  data. 

The  curves  in  Fig.  1  also  suggest  that  the  ideal  window  size  for  both  algorithms  is  in  the  upper  teens  or  low  twenties. 
Though  some  of  the  highest  detections  were  above  this  range,  they  were  not  dependably;  the  total  number  of 
detections  actually  declined  as  the  window  size  approached  49x49  for  three  and  four  of  the  nights  with  the 
Gaussian-PSF  and  Gaussian-Dirac  algorithms,  respectively.  Also,  neglecting  efficiencies  that  can  be  found  with 
vector  processing,  the  computational  cost  of  computing  either  algorithm  increases  roughly  proportional  to  the  area 

of  the  window,  N 2 . 

As  can  be  seen  in  fig.  2,  the  Gaussian-Dirac  algorithm  does  register  an  appreciable  amount  of  detections  that  the 
Gaussian-PSF  algorithm  does  not.  If  forced  to  choose  only  one,  the  Gaussian-PSF  might  yield  better  results; 
however,  the  best  option  might  be  to  choose  both.  If  the  same  window  size  is  used  for  both  algorithms,  then  the 
computational  burden  is  significantly  lower  since  the  window  is  already  in  memory  and  the  mean  has  been 
computed. 


8.  CONCLUSION 

This  paper  has  compared  the  Gaussian-Dirac  and  the  Gaussian-PSF  log-likelihood  detection  algorithms  and  shown 
how  search  grid  size  and  window  size  affect  their  performance  in  detecting  a  satellite  in  varying  degrees  of 
illumination.  The  Gaussian-PSF  performed  better,  as  was  expected.  This  study  should  help  in  selecting  the  optimal 
parameters,  whichever  algorithm  is  used  in  the  detection  of  dim  objects. 
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